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1. Introduction 

The calculation of Ruelle-Pollicott resonances is a nontrivial problem in dynamical 
systems theory from both a theoretical as well as a numerical viewpoint. They describe 
the asymptotic decay of phase space densities towards the invariant phase space measure. 
The Perron-Frobenius (PF) operator gives the propagation of phase space densities, and 
the Ruelle-Pollicott resonances are the eigenvalues of the coarse-grained PF operator. 
The coarse-graining is a crucial element in their description. In a chaotic system, the 
stretching and folding of phase space structures causes the phase space densities to 
attain finer and finer structures under the evolution of a PF operator. It is only after 
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the coarse graining that these structures are smeared out and that the decay towards 
the invariant measure can be observed. 

Here we are concerned with numerical approaches towards the calculation of Ruelle- 
Pollicott resonances. There exist two typical approaches towards their calculation. The 
most common one is to express the PF operator in a finite basis that also provides a 
natural way of coarse graining, see e.g. [121 H3]- This reduces the problem to the one 
of finding the eigenvalues of a matrix. The problem with such an approach is that it 
requires large matrices in order for the resolution to be sufficient for good convergence. 
With increasing matrix size, however, the coarse graining is also reduced, causing the 
creation of a large number of spurious eigenvalues that, in the case of area preserving 
invertible maps, tend to fill out the whole unit circle. This makes even the identification 
of resonances a difficult task [3]. This can be considered as a many-observable, single- 
timestep approach. The numerical difficulties stem from the fact that, in order to 
obtain a good long-time description of the system, a large basis needs to be considered 
to faithfully capture the actual dynamics from a single timestep. 

A different approach considers taking only a single observable (phase space function) 
and, instead of focusing on a single iteration of the PF operator, the auto-correlation 
function of such an observable is computed over many steps of the iteration. By using 
linear predictors or Pade approximates the signal is fitted to a sum of exponentials, 
which in turn give the positions of the resonances, see e.g. !8j . A similar technique is 
used in jHJ [15] for determining the spectral lines of quantum systems. A whole class 
of numerical schemes based on this approach is presented in [TT] . Compared to matrix 
diagonalization, this approach is computationally cheaper. The disadvantage of such an 
approach is a certain amount of arbitrariness when choosing the number of exponentials 
to be fitted to the signal. As will be shown, this can introduce significant uncertainties 
into the obtained results. Also, while such an approach faithfully captures the long time 
behaviour of the system, the results also depend on the behaviour of a single chosen 
observable. This might cause certain resonances to be missed or poorly resolved. 

We present a hybrid method that combines the favourable aspects of the described 
approaches. It uses a moderate number of observables as well as a moderate number 
of timesteps when estimating the eigenvalues of the PF operator. It will be shown that 
such an approach requires fewer observables than the direct diagonalization technique, 
and it provides its own criterion on the optimal number of time steps. By using multiple 
simultaneous observables, it also avoids the problem of missing resonances. 

In the following, we define the Perron-Frobenius operator and the Ruelle-Pollicott 
resonances. We then discuss the issues with the linear predictor approach and give 
an analytical estimate on the spurious "ring" of resonances that can cause numerical 
difficulties. The hybrid method is presented next, followed by the comparison of results 
for various methods, and conclusions. 
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2. Preliminaries 

We consider an ergodic and mixing [TJ time-discrete dynamical system /* : S — > S in 
the phase-space S with the invariant measure /i : S — > IR + normalized so that p(S) = 1. 
A single time step of the system is defined by the iteration 

x t +i = f{x t ) x t e S . (1) 

Instead tracking individual trajectories the dynamics can be expressed by the evolution 
of the their probability distribution over the phase space p l : S — > R+ defined using the 
Perron-Frobenius (PF) operator L as 

p t+1 (x) = (Lp t ){x) = J dy5(x - f(y)) Pt (x) . (2) 

This can can be explicitly written in the form of a sum over all y that are mapped into 
x by the map /(•) reading 

yef-Hx) 1 w vy;;i 

where f'(x) is the Jacobian of the map. For area preserving invertible maps the invariant 
measure is Lebesgue dfi(x) — dx/ J s dy and the action of the PF operator reduces to 

p t+1 (x) = (L Pt )(x)= Pt (r 1 (x)). (4) 

Let us denote by p(g) the phase-space average J dp(x)g(x) of a function g(x) w.r.t. to 
the invariant measure p(-). Then the correlation of two real smooth observables u(x) 
and v(x), with their phase-space averages p{u) and p{v) equal to 0, is given by 

C u , v (t) = p{u ■ L l o v) . (5) 

In mixing systems and under certain conditions [2J [3] , it has a well defined asymptotic 
time dependence of the form 

c UiV (t) ~ afHA^*/} , (6) 

where A u ^ v is a constant depending on used observables and v is called the leading 
Ruelle-Pollicott (RP) resonance. The RP resonances are the eigenvalues of the coarse 
grained PF operator [3]. Here we are interested mostly in dynamical systems confined 
to a compact phase space in which the RP resonances are believed to depend on the 
fine structure of the mixing process as opposed to the diffusion in phase space. Namely, 
for a large area preserving maps on the cylinder (an infinite phase-space) it was shown 
in [HI [Tj that there exists a strict connection between the resonances and the diffusion 
process. 



3. Discussion of the linear predictor approach to Ruelle-Pollicott 
resonances 

A common type of approach towards calculating the resonances of arbitrary signals is 
to use linear predictors, which match the signal (or rather its auto-correlation function) 



A hybrid method for calculation of Ruelle-F 'ollicott resonances 



4 



to a sum of exponentials [El [111 US, [UJ. If we consider a dynamical system with 
the PF operator L, then in the linear prediction (LP) approach we assume that an 
auto-correlation function C(t) of some observable u(x) with a zero phase-space average 
fi{u) = can, for t > 0, be represented by a finite sum of exponentials as 

r 

C u (t)= f i(u-L t u)^J2 C ^^ ( ? ) 

i=l 

where {^}[ =1 are the LP resonances and Cj 6 C are constant weights associated with 
the individual resonances. In an ideal case, the LP resonances represent or approximate 
the RP resonances of the dynamical system. Equivalently, the auto-correlation signal 
C(t) can approximated by a finite linear predictor (LP) [TU IT?] , which is given in the 
form of a finite difference recursion 

r 

C(t) = ^d i C(t-i)+^t), t = r,r + l,... , (8) 

i=l 

with a discrepancy £ t between the value produced by the finite LP model and the true 
value of the auto-correlation function C(t). There are several known methods to obtain 
LP coefficients di [IB], from which the least-squares is the most straightforward and 
therefore used here. 

For a given auto-correlation signal {C(t)}" =0 of length n, the vector of LP 
coefficients d = (di) r i=1 is, in the least-squares approach, chosen such that the total 
discrepancy X]™=t-£(^) 2 * s minimal. This condition yields a linear system of equations 
for the vector d written as 

Ad = b, (9) 

with the matrix A = [Aj]ij=i an d the vector b = \bj\ r i=l given by 

n—r n—r 

A iJ = y ^2c(t + r-i)C(t + r-j), k = ^C(t + r - i)C{t + r) . (10) 
t=o t=o 
This system of equation is usually well conditioned and the solution - the vector d - 
determines the form of the LP. The LP resonances {^j}[ =1 are given as the roots of the 
characteristic polynomial 

r 

P LP (x) = x r -Y,d i x r - 1 (11) 
i=i 

corresponding to the LP. 

The number of exponential terms r chosen is a critical parameter in such 
approaches. By choosing too few exponentials, the positions of the resonances are 
not properly resolved. Here we will show that, on the other hand, by taking too many 
exponentials, noise causes the fitting procedure to try and follow a random signal, which 
pushes the calculated resonances out towards the unit circle. We demonstrate this by 
considering one realization of a typical auto-correlation function of a real observable in 
system with only two clearly distinguishable complex resonances v and v* written as 

C(t) = K{z/} + e(t) , t = 0,1,2,.... (12) 
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The function e(t) represents a noise resulting from a finite precision of the auto- 
correlation measurement or calculations. In the case r > 2 and e(t) = the auto- 
correlation possesses only two non-zero LP resonances v and v* . The LP resonances 



corresponding to one realization of the auto-correlation (12) with noise is depicted in 
figure [1} We may see that the true resonances - the complex conjugated pair v and v* 
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Figure 1. The spectrum of LP resonances corresponding to the the auto-correlation 
C(t) = 5i{z/} + e(i) for a different number of considered resonances r with v = 0.3e l0 ' 5 , 
t E [0,n = 10 3 ] and e(t) being random and uniformly distributed in [— 10~ 12 , 10 -12 ]. 
The radius of the circle going through the ring of resonances has been estimated by 
R(r — 2, n) (17). The gray-shaded discs mark the positions of v and v* . 



are present in the set of LP resonances. Beside the pair, there is a ring of resonances 
which is moving towards the unit circle as the number of resonances r is increased. The 
ring of resonances is present whenever the number of considered resonances r exceeds 
the number of detectable resonances. In the following we will present a model giving a 
qualitatively correct radius of the ring of LP resonances. 

Considering a LP model using r resonances, we see that the ring radius composed of 
q resonances is virtually independent from the amplitude of r — q detectable resonances. 
We expect this ring to occur due to the noise present in the calculation of the correlation 
functions. This radius can therefore be expected to be estimated from the LP resonances 
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corresponding to a numerical noise auto-correlation model C mo dei(£) = <5t,o + for 
t G [0, n] using q resonances, where we assume that e(t) is random white noise with 
variance v = (e(t) 2 ) much smaller than 1. By inserting C mo d e \(t) into equations (10) 
and taking into account the central limit theorem |9] , we find that the matrix elements 
Ai j and the vector elements 6j are, in the limit fi>l, approximately Gaussian variables 
with statistical moments 

(Aij) = 5 itj [5 hq + v{n - q + 1)] , (bi) = 0, {bj) = v[6 i>q + v(n - q + 1)] , (13) 

where i,j G [1, q] and the brackets (•) represent the averaging over different noise 
realizations. The average of the matrix A is diagonal and by neglecting the off-diagonal 
elements one can approximate the LP coefficients by 

dtM-j-, (14) 

which means that di are approximately Gaussian variables with a zero average (di) = 
and variances equal to 

<<*?> = S ^ T\ ■ ( 15 ) 

x l/ 8 i>q + v{n - q + 1) 

In the limit of small noise e(t) — > 0, the coefficient d q = 0, and the remaining non-zero 
LP coefficients (d\, . . . , d q -i) define a characteristic random polynomial 



LP, model 



X 9 " 1 - 2^ diX q - l j (16) 

corresponding to different realizations of the auto-correlation C mo d e i(t)- 

We find empirically that the roots of the polynomials x s — £ YltZo a i x% °f or der s, 
with coefficients a, being normalized Gaussian variables iV(0, 1) and |^| 1, form rings 
in the complex plane with the radius approximately £ » . This radius can also be obtained 
if we consider the zeros of a polynomial containing only the highest and constant order 
terms of the original polynomial, x s — £ = 0, with the constant term replaced by its 
root mean squared value. According to this rather ad-hoc expression, the radius of 
resonances corresponding to C mo d c i(t) is approximately given as 

R(q,n) « (n - q+ 1)"^^) . (17) 

The comparison of this estimate with the average amplitude of the largest LP resonance 
is shown in figure [2} We see that the analytic estimate is not very accurate, but it gives 
a qualitatively correct functional dependence. Similar estimates for the average radius 
of the ring of resonances in the complex plane can also be found for other approaches to 
obtain the LP coefficients, such as the auto-correlation method (sometimes called the 
Yule- Walker equation), the Burg's method etc. 

While having the benefit of detecting the asymptotic behaviour directly, the linear 
predictor method is therefore disadvantaged by spurious solutions that make it hard to 
distinguish the real resonances from noise in practice. Furthermore, by focusing on a 
single observable, the results may significantly depend on the choice of the observable, 
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Figure 2. The average amplitude of the largest LP resonance of the model auto- 
correlation C(t) — So,t + e(i) for different number of considered resonances q and 
lengths n of the auto-correlation signal. The inset shows the difference A between the 



average amplitude of the model and their analytic estimate (17). The e(t) is random 
and uniformly distributed [— 1CT 12 , 10~ 12 ]. 



which is nicely demonstrated in [10]. On the other hand the LP methods are fast. The 
time complexity of preparing the data for the LP methods is 0(n), while calculating 
the LP resonances in a robust way takes 0(r 3 ) times steps. 

4. The hybrid method of detecting Ruelle-Pollicott resonances 

Here we propose a hybrid method that relies on the diagonalization of relatively small 
matrices, yet takes into account the information of multiple iterations of the PF operator. 
A somewhat analogous approach was taken in [18] in the analytic exploration of the 
resonance structure of the kicked rotor. 

We start with the same approach as in the case of the direct diagonalization, but 
instead of a single matrix we compute a time-sequence of matrices. Let us choose a basis 
B r = {(fii}l =l for the densities, which is orthonormal w.r.t. to the invariant measure 
n((pi ■ (fj) = The elements of the matrix corresponding to the n-th iteration of the 
map is given by 

if = n (cp* ■ L* o Vj ) . (18) 

One may also interpret this as a matrix of various cross correlations. In the usual 
diagonalization approach, one takes the matrix and computes its eigenvalues to 
obtain the resonances. In the linear prediction and related approaches, the diagonal 
matrix elements are taken and fitted with a sum of exponentials as a function of t. 

Our approach is based on the principle that the iterations of the exact PF operator 
are related by 

L t+1 = L o V . (19) 
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The subsequent matrices T® should therefore be related by a single t-independent 
matrix Y in an approximate way, 

T (t+i) w Y T® _ (20) 

The goal is to find a procedure to calculate the matrix Y such that the above relation 
is best satisfied for all t. The eigenvalues of such a matrix are then expected to give the 
RP resonances. A way to calculate an optimal Y is to minimize the error norm 

n-1 

E(Y,Y^) = Y^WtTr [(T®W - T^) (YT® - T^)] . (21) 
t=o 

The weights are introduced in order to boost the information from the later time steps, 
which provide the most information about the asymptotic decay, but would otherwise 
be drowned by the much larger correlations found in the first few timesteps. The choice 
is 

w t = [Tr (T^T®) - T ]~\ (22) 

where 7 = if the basis set does not includes the uniform density and 7 = 1 if the basis 
set can fully reproduce the uniform density, as this is a trivial contribution that can 
be subtracted. With such a choice, information from all the time steps is treated with 
roughly the same significance. 

We may treat the matrix elements of the matrix Y and its adjoint as independent 
complex numbers. We then try and minimize E with respect to all the matrix elements 
by setting its gradient to 0. One may note that 

' Tr(AB) = —— V" V" A kl B lk = Bji . (23) 



The requirement that the gradient of Tr AB with respect to all elements of A equals 
then reduces to the matrix equation B = 0. As the trace operation is invariant with 
respect to cyclic permutations of matrix products, we rearrange the order of matrix 



multiplications in equation (21) to obtain 

n-1 

E(Y, Ft) = Y,w t Ti [r f (YT®T®^ - T^+VtW) + T {t+l ^ (YT® - T^)] (24) 
t=o 

and equaling the gradient of g with respect to the elements of Y^ to then leads to the 
equation 

n-1 

w t (YT®T<M - T^T^t) = . (25) 

t=o 

Solving this for Y yields 

Y= (j2^tT (t+1) T®^ ^w t T®T®^j . (26) 
The RP resonance can be identified as the (sub-) dominant eigenvalue of the matrix Y. 



A hybrid method for calculation of Ruelle-F 'ollicott resonances 



9 



It should be noted here that the presented linear predictor class methods can 
be extended so that not only a single auto-correlation function signal is fitted by 
exponentials, but signals of cross correlations of a set of observables given by equation 



(18) are considered. In such an approach, if the predictor is only a single step one, a 
super-matrix Q is sought such that the relationship 

k,l 

is best satisfied for all iterations. If matrices are considered as vectors, the 
eigenvalues of Q are then expected to give the resonances of the system. With r giving 
the basis size of the system, we obtain r 2 eigenvalues, whereas the direct diagonalization 
and hybrid approaches towards the calculation of the RP resonances return r eigenvalues. 
This suggests a certain amount of underdetermination of the problem. It is likely for this 



reason that, when tested, the approach as given in equation (27) proved very unstable 
with a lot of spurious eigenvalues. 



Equation (20) can be considered as a particular reduction of equation (27) by only 



taking a very specific connection between correlations from subsequent timesteps. This 



reduction reflects the actual underlying process as given by equation ( 19 ) and eliminates 
the underdetermination of the problem. 

The calculation of all the matrices {T^}™ =0 to some precision is a lengthy process 
having the time complexity of at least 0(nr 2 ), but typically even more (0(nr 3 ) for the 
two-dimensional example studied later) as the integration grid has to be refined for faster 
oscillating observables, whereas obtaining the effective propagator Y from matrices T" 
and calculating its eigenvalue spectrum takes 0(r 3 ) steps. In all usual cases the time 
complexity of direct diagonalization (n = 1) and of the hybrid method n > 1 is several 
orders larger than for the LP methods when using the same number of resonances and 
time steps considered. 

5. Results of the hybrid method 

The presented method for searching of the RP resonances was testing using an area 
preserving map - the perturbed Arnold's cat map defined on a two-dimensional torus 
T 2 = [0, l) 2 (19] with one iteration step given by 

K 

y t+l = y* + x l - — sin(27rx*) (mod 1) , (28) 

2n 

x t+1 = x l + y t+1 (mod 1) . (29) 

This map can be shown to be fully ergodic for K £ [0, 1]. At K — 1 it develops a 
marginally stable periodic orbit at the origin. We therefore expect some RP resonances 
of the system to approach unity in terms of their absolute value when approaching 
K — > 1. We choose to work with the real Fourier basis defined as 



1 



y?n(^) = { \/2cos(27mx) 
a/2 sin(27rna;) 



n = 

n = even > (30) 
n = odd 
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because they are smooth, point convergent and non-local. We use p Fourier modes along 
each axis and so the functions on the torus are spanned by the functional basis 

£ P 2 = {fiA x > v) = <Pi(x)vj(y)} p ijlo ■ ( 31 ) 

Notice that the dimension of the vector space that we work with is r = p 2 . 




Figure 3. The dependence of the leading RP resonance for the Perturbed cat map 
system at K = 0.9 on the trajectory length n and the number of base functions spanned 
along one axis p, using the hybrid method with 7=1 (a)-(c), and the evolution of 
the trace Tr(T^T^ for p = 20 (d). The overlap integrals represented by the matrix 
elements T^j are calculated by extrapolating their Simpson's rule approximants |20j 
on 1000 x 1000, 2000 x 2000 and 4000 x 4000 grids of points. 

In figure [3] we show the results of our approach for the case of K = 0.9, where the 
system is already close to being marginally stable. We show the absolute value of the 
sub-dominant eigenvalue as a function of the basis size p and the number of timesteps n 
taken. The case of n — 1 is equivalent to a direct diagonalization approach, whereas the 
case of small p corresponds more closely (but is not equivalent) to the linear predictor 
methods. As we can see in figure [3]^a), if we only increase p while keeping n = 1 
(standard diagonalization), or if we keep p = 1 and increase n, the convergence is poor. 
Good convergence is only obtained by increasing both p and n simultaneously, which 
is also the main idea behind the hybrid approach. Figures |3](b) and [3|c) demonstrate 
the convergence properties when taking either a constant p (b) or constant n (c) cross 
section. It is interesting to note that direct diagonalisation (n = 1) is very sensitive 
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Figure 4. The leading RP resonance as a function of auto-correlation signal 
C(t) = fi(g ■ L l o g) of length n and considered resonances r obtained using the linear 
predictor from the auto-correlation of an observable g(x,y) = A[sin(7rx 2 ) cos(27ry) — 
sin(27rx) cos(27ry 2 )] in the perturbed cat map system for K = 0.9 with A chosen so that 
Jj x » 2 dxdy g(x,y) 2 = 1. The C(t) is calculated with the absolute precision 3 • 10~ 6 

by averaging the time auto-correlations C(i,£) = X)"=o* _1 ff(0s(/~*(£)) °f length 
n = 1024 over 10 8 initial points uniformly distributed on the torus £ € [0, l) 2 . 



to the parity of the basis size p, but that the hybrid method tends to even out these 
differences with increasing n. 

Figure |3^b) also shows a significant change in convergence for roughly n > 10. In 
figure ^d), the weight as given in equation (22) is shown. In an exact calculation this is 
expected to drop to 0, however, due to the noise present in the correlation computation, 
the number saturates at a higher level. The time n at which the saturation is reached 
is also the cutoff point for the hybrid calculation, as is reflected in the convergence 
properties of figure [3](b) . 

In figure [4] we give the largest resonance as obtained via the linear predictor method 
for an observable chosen to not have any of the system symmetries. We can see that 
the procedure converges when increasing the number of chosen resonances r, until a 
certain value is reached, beyond which its value tends towards the unit circle. This is 
representative of the phenomenon as explained in section [3] and is a serious limitation 
of the single variable linear prediction methods. 

In figure [5] we show a comparison of full resonance spectra obtained via linear 
predictors in the least square approach and the hybrid method proposed here. We see 
that using the direct diagonalization (n = 1) and the hybrid method with n — 10 
they give essentially different spectra for same number of spectral points given p 2 , with 
the direct diagonalization tending to an underestimation of the main resonances. The 
spectrum obtained via LP captures the main resonance well and agrees with the hybrid 
method, however even for r = 10 resonances one can observe that it is only the main 
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Figure 5. The spectrum of resonances obtained using the hybrid method and the LP 
approach for different time intervals [0, n] and the number of considered resonances r 
in LP, or basis size along one axis p in the hybrid method for the perturbed cat map 
system with K = 0.9. For the details on numerical calculation of matrix elements T^j 
and used auto-correlations see the captions under figures [4] and [3j 

resonance that is attainable by the LP approach, as all the other calculated resonances 
that roughly form a circle are representative of the phenomenon described in section [3] 

The convergence properties of the hybrid method for the full spectrum as a function 
of the basis size are demonstrated in figure [6] It can be seen that for all values of 
p = 7, 10, 15, the leading eigenvalues already appear to have converged quite well for 
n = 10, this being the natural cutoff time as explained earlier. For all these cases, there 
also exists a cloud of eigenvalues around the origin with a roughly constant radius that 
likely corresponds to the essential spectrum of the PF operator. Its radius slowly grows 
with p, which is believed to be due to the noise in the matrix element computations, a 
phenomenon not entirely unlike the one observed in the LP approach, but which does not 
affect the calculations in a critical manner. For comparison, the direct diagonalization 
results are also given. 

As an interesting example, we may consider the case of the perturbed cat map for 
K = 1.1. At this parameter, the map ceases to be ergodic and attains regular islands. 
The " resonances" corresponding to the regular islands should have their absolute values 
at unity Figure [7] repeats the calculation of figure [3] for such a system. We may see 
that the hybrid method easily captures these "resonances" close to unity, whereas both 
direct diagonalization as well as linear predictors (not shown) struggle to reach this 
value and underestimate the "resonances" . One interesting failure of the method can be 
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Figure 6. The spectrum of resonances obtained using the hybrid method for different 
basis sizes along one axis p for the perturbed cat map system with K = 0.9. The 
direct diagonalization approach (n — 1) for p — 7 is shown as crosses for comparison. 
For the details on numerical calculation of matrix elements T^} and auto-correlations 
used see the captions under figure [3| 



seen for large p at around n — 10, where the largest eigenvalue is overestimated and is 
pushed outside the unit circle. As can be seen in figure ff\d), the weight w n never drops 
towards but maintains a plateau after n > 10 due to the regular islands. It is only 
after this time that the initially larger chaotic signal decays and the small regular islands 
are then properly resolved in the correlation functions. Around n m 10 there exists a 
competition between the chaotic signal contributions and the regular ones, and this is 
the likely reason why the hybrid method fails for when terminating the calculation at 
precisely these n. 



6. Conclusions 



We have presented a new numerical method to calculate the RP resonances in dynamical 
systems with compact phase-space, which is particularly useful when one is limited to 
a small number of observables. The method is a hybrid between the linear prediction 
and diagonalization of the coarse-grained Perron-Frobenius operator. It considers many 
simultaneous observables, which is the main advantage of the direct diagonalization 
technique, over many time steps, which is the advantage of the linear prediction 
methods. It also requires far less observables for a good convergence than the direct 
diagonalization, and provides its own criterion as to the number of time steps that give 
best results. It is shown in the perturbed cat map system that the presented method, in 
the contrast to other mentioned approaches, gives a convergent and unambiguous value 
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1.25 




Figure 7. The leading RP resonance of the perturbed cat map at K = 1.1 as a 
function of the trajectory length n and the number of base functions spanned along 
one axis p using the hybrid method with 7=1 (a)-(c) and the evolution of the trace 
Tr(T^) 2 ) for p = 20. Numerical details of calculations are in the caption under figure 

13 



for the RP resonances. 
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